Reexamination of Hagen— Poiseuille flow: 
shape-dependence of the hydrauHc resistance in microchannels 
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We consider pressure-driven, steady state Poiseuille flow in straight channels with various cross- 
sectional shapes: elliptic, rectangular, triangular, and harmonic-perturbed circles. A given shape is 
characterized by its perimeter V and area A which are combined into the dimensionless compactness 
number C — V'^/A, while the hydraulic resistance is characterized by the well-known dimensionless 
geometrical correction factor a. We find that a depends linearly on C, which points out C as a 
single dimensionless measure characterizing flow properties as well as the strength and effectiveness 
of surface-related phenomena central to lab-on-a-chip applications. This measure also provides a 
simple way to evaluate the hydraulic resistance for the various shapes. 

PACS numbers: 47.60.-|-i, 47.10.-|-g 
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I. INTRODUCTION 

The rapid development in the field of lab-on-a-chip 
systems during the past decade has put emphasis on 
studies of shape-dependence in microfluidic channels. 
Traditionally, capillary tubes would have circular cross- 
sections, but today microfabricated channels have a va- 
riety of shapes depending on the fabrication technique 
in use. Examples are rectangular channels obtained by 
hot embossing in polymer wafers, semi-circular channels 
in isotropically etched surfaces, triangular channels in 
KOH-etched silicon crystals, Gaussian-shaped channels 
in laser-ablated polymer films, and elliptic channels in 
stretched PDMS devices, see e.g., Ref.hl 

The pressure-driven, steady-state flow of a liquid 
through long, straight, and rigid channels of any constant 
cross-sectional shape is referred to as Hagen-Poiseuille 
(or simply Poiseuille) flow, and it is often characterized 
by the hydraulic resistance, -R^yd = ^p/Q^ where Ap 
is the pressure drop along the channel and Q the flow 
rate through the channel. In Fig. ^ is shown an arbitrar- 
ily shaped cross-section in the xy plane for a straight 
channel placed along the z axis. A natural unit for the 
hydraulic resistance is given by dimensional analysis as 




FIG. 1: An arbitrary cross-sectional shape Q with perimeter 
dQ of a straight fluid channel with pressure-driven steady- 
state flow. The contours show the velocity v{x,y) obtained 
numerically from Eq. (|^ by a finite-element method. The 
velocity is zero at the boundary and maximal near the centre- 
of-mass. 



R^y^ = rjL/A^, where L is the channel length, rj the dy- 
namic viscosity of the liquid, and A = Jj^ dxdy the cross- 
sectional area. Typically, the fluid flow is subject to a 
no-slip boundary condition at the walls 9f2 and thus the 
actual hydraulic resistance will depend on the perimeter 
as well as the cross-section area. This dependence can 
therefore be characterized by the dimensionless geomet- 
rical correction factor a given by 

a=T^. (1) 



^hyd 

In lab-on-a-chip applications 0, Q , where large surface- 
to-volume ratios are encountered, the problem of the 
bulk Poiseuille flow is typically accompanied by other 
surface-related physical or bio-chemical phenomena in 
the fluid. The list of examples includes surface chemistry, 
DNA hybridization on fixed targets, catalysis, interfacial 
electrokinetic phenomena such as electro-osmosis, elec- 
trophoresis and electro-viscous effects as well as contin- 
uous edge-source diffusion. Though the phenomena are 
of very different nature, they have at least one thing in 
common; they are all to some degree surface phenomena 
and their strength and effectiveness depends strongly on 
the surface-to-volume ratio. It is common to quantify 
this by the dimensionless compactness C given by 



(2) 



where V = Jg^ d£ is the perimeter of the boundary 9f2 
confining the fluid, see Fig. ^ For other measures of C 
we refer to Ref. Q and references therein. In this paper 
we demonstrate a simple dependence of the geometrical 
correction factor a on the compactness C and our results 
thus point out a unified dimensionless measure of flow 
properties as well as the strength and effectiveness of 
surface-related phenomena central to lab-on-a-chip ap- 
plications. Furthermore, our results allow for an easy 
evaluation of the hydraulic resistance for elliptical, rect- 
angular, and triangular cross-sections with the geometri- 
cal measure C being the only input parameter. Above we 
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have emphasized microfluidic flows because here a vari- 
ety of shapes are frequently encountered. However, our 
results are generally valid for all laminar flows. 



II. POISEUILLE FLOW 

Due to translation invariance along the z axis the ve- 
locity field of a Newtonian fluid in a straight channel is 
parallel to the z axis, and takes the form v = v{x^y)e2.- 
Consequently, the non-linear term in the Navier-Stokes 
equation drops out 4], and in steady-state, given the 
pressure gradient — (Ap/i)ez, the velocity v{x, y) is thus 
given by the Poisson equation. 



Ap 



(3) 



with the velocity being subject to a no-slip condition at 
the boundary dil. The relation between the pressure 
drop Ap, the velocity v{x,y), and the geometrical cor- 
rection factor a becomes 

Ap = RhydQ = a-Rhyd<3 = "Kyd dxdy v{x, y), (4) 

where Q is the volume flow rate. 
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FIG. 2: Correction factor versus compactness for the ellipti- 
cal, rectangular, and triangular classes. The solid lines are 
the exact results, and the dashed lines indicate Eqs. Q, l)14|l. 
and list . Numerical results from a finite-element simulation 
are also included (o, □, and A). Note that in the case of 
triangles all classes (right, isosceles, and acute/obtuse scalene 
triangles — marked by different grayscale triangles) fall on 
the same straight line. 



III. THE GEOMETRICAL CORRECTION 
FACTOR VERSUS COMPACTNESS 

Our main objective is to find the relation between the 
geometrical correction factor a and the compactness C 
for various families of geometries. 



A. Elliptical cross section 

The elliptical family of cross-sections is special in the 
sense that Eq. ^ can solved analytically (see e.g. Ref.|J) 
and we can get an explicit expression for the geometrical 
correction factor introduced in Eq. (^). For an ellipse 
centered at the origin with semi-major and minor axes a 
and b it can be verified by direct insertion that 



v{x,y) 



Ap {abf 
^2(a2-f 62) 



fulfils Eq. From Eq. Q it can now be shown that 

a(7) = 47r(7 + 7-i) (6) 
where 7 = a/6. Furthermore, for an ellipse we have 



(7) 



The relation between a and C can now be investigated 
through a parametric plot. In order to get an approxi- 
mate expression for a(C) we begin by inverting Eq. I^. 



By selecting the proper root we get 7(a) which we then 
substitute into Eq. ([7|l such that 



a2 - (87r)2cos6' 



(8) 



Expanding around a = Stt and inverting we get 

a(C) = ^C-|^ + 0([C-4^P), (9) 

and in Fig. [21 we compare the exact solution (solid line), 
from a parametric plot of Eqs. © and 0, to the approx- 
imate result (dashed line) in Eq. (j^ . Results of a numer- 
ical finite-element solution of Eq. © are also included (o 
points). As seen, there is a close-to-linear dependence of 
a on C as described by Eq. @ . 



B. Rectangular cross section 

For a rectangle with width-to-height ratio j = w/h we 
solve Eq. Q using Fourier series jj]. 



ApAh^ 
v[x,y) = —r — 
•qL TT"^ 



5Z n3 ( "'^ 



n=l,3,5,. 



cosh(ri7rx//i) 
\ cosh(n7ru'/2/i) 



(10) 



sin(n7rj///i) 



is indeed a solution. Here, the coordinate system is cho- 
sen so that —w/2 < X < w/2 and < y < h. From 



3 



Eq. Q it follows that 



«(7) = ^ 

\n=l,3,5,. 



717 



■ tanh(n7r7/2) 



and for the compactness we have 

C(7) =8 + 47 + 4/7. 
Using that tanh(x) ~ 1 for x 3> 1 we get 
127r572 



a(7) 



7r57 - 186C(5) ' 



7> 1, 



(11) 
(12) 

(13) 



and substituting 7(C) into this expression and expanding 
around C(7 = 2) = 18 we get 



a{C) 



22 



■C 



65 

y 



0([C - 18] = 



(14) 



For the two Taylor coefficients we have used the first three 
terms in the continued fraction. In Fig.|21we compare the 
exact solution, obtained by a parametric plot of Eqs. 
and (|12|1 . to the approximate result, Eq. (|14|l . Results of 
a numerical finite-element solution of Eq. JSJ are also 
included (□ points). As in the elliptical case, there is 
a close-to-linear dependence of a on C as described by 
Eq. 



C. Triangular shape 

For the equilateral triangle it can be shown analytically 
that a = 20^/3 and C = 12^/3, see e.g. Ref. 0. However, 
in the general case of a triangle with side lengths a, b, 
and c we are referred to numerical solutions of Eq. 
In Fig. 121 we show numerical results (A points), from 
finite-element simulations, for scaling of right triangles, 
isosceles triangles, and acute/obtuse scalene triangles (for 
the definitions we refer to Ref. la). The dashed line shows 
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(15) 



where the slope is obtained from a numerical fit and sub- 
sequent use of the first three terms in the continued frac- 
tion of this value. As seen, the results for different classes 
of triangles fall onto the same straight line. Since we have 



C(a, 5, c) = 



8{a + b + c)' 



^i(a2 + 62_,_c2)'- (a4 + 64 + c4) 



(16) 



the result in Eq. I|15|l allows for an easy evaluation of 
-Rhyd for triangular channels. 




(b) 



X = p cos { 
y = p sin 6 



r = ap[l + esin(fc(?)] 

X = ap\l + e sin(fc6)] cos ( 

y = ap^l + esin(A;9)] sin 6 



FIG. 3: (a) The geometry of the unperturbed and analyti- 
cally solvable cross section, the unit circle, described by co- 
ordinates {S:,y) or {p,9). (b) The geometry of the perturbed 
cross section described by coordinates {x,y) or (r, 0) and the 
perturbation parameter e. Here a = 1, k = 5 and e = 0.2. 



D. Harmonically perturbed circle 

By use of shape perturbation theory it is possible to 
extend the analytical results for Poiseuille flow beyond 
the few cases of regular geometries that we have treated 
above. In shape perturbation theory the starting point 
is an analytically solvable case, which then is deformed 
slightly characterized by some small perturbation param- 
eter e. As illustrated in Fig. |3|the unperturbed shape is 
described by parametric coordinates (5;,y) in Cartesian 
form or (p, 9) in polar form. The coordinates of the phys- 
ical problem we would like to solve are (x, y) in Cartesian 
form and (r, cj)) in polar form. 

As a concrete example we take the harmonic pertur- 
bation of the circle defined by the transformation 



r = a p [1 + e sin(A:6')] , 
x(p, 6) = a p[l + e sm{kd)\ cos 0, 
y{p, 6*) = ap[l + esin(A;6')] sin 6*, 



(17a) 
(17b) 
(17c) 
(17d) 



where a is length scale, k is an integer (> 2) defining 
the order of the harmonic perturbation, < < 27r, and 
< p < 1. For e = the shape is unperturbed. The 
boundary of the perturbed shape is described by fixing 
the unperturbed coordinate p = 1 and sweeping in 9, 



an 



(x,y 



ihO], y[l,9]). 



(18) 



It is desirable to formulate the perturbed Poiseuille prob- 
lem using the unperturbed coordinates. To obtain ana- 
lytical results it is important to make the appearance of 
the perturbation parameter explicit. When performing a 
perturbation calculation to order m all terms containing 
e' with I > m are discarded, while the remaining terms 
containing the same power of e are grouped together, and 
the equations are solved power by power. To carry out 
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the calculation the velocity v{x,y) is written as 

v{x,y)^v{x[p,e],y[p,e]) (19) 



(p, e) + ev^'>{p, 0)+e' v^'>{p,9) + --- 

Likewise, the Laplacian operator in Eq. must be ex- 
pressed in terms of p, 9, and e. The starting point of this 
transformation is the transformation of the gradients 



dr = {drp) dp - 
d<t, = {d4,p) dp 



(20a) 
(20b) 



The derivatives {dj.p), (d^O), {d^p), and {d^9) are ob- 
tained from the inverse transformation of Eqs. (|17a|l and 
(|17b|l . The expansion in Eq. (|19|l can now be inserted 
into Eq. (0) and using the derivatives, Eqs. H20a() and 
(|20b|l . we can carry out the perturbation scheme. The 
calculation of the velocity field to fourth order is straight- 
forward, but tedious. With the velocity field at hand we 
can calculate the flow rate and from Eq. Q we get 



l + 2(fc-l)e2 



47 - 78fc + 36fc2 - 4fc3 
8 



(21) 



where we have used the exact result A — + ) ^o, 
for the area. The result only involves even powers of e 
since e — > — e is equivalent to a shape-rotation, which 
should leave a invariant. From an exact calculation of 
the perimeter V we get the following expression for C, 



C = 47r + 27r(fc2 - 1) 



(22) 



Since a is also quadratic in e this means that a depends 
linearly on C to fourth order in e, 



a(C) 



1 + fc 



C-8^.-HO(.^). (23) 



Note that although derived for fc > 2 this expression 
coincides with that of the ellipse, Eq. ||5J), for k = 2. 
Comparing Eq. (|21|l [to second order in e] with exact 
numerics we find that for e up to 0.4 the relative error is 
less than 0.2% and 0.5% for fc = 2 and fc = 3, respectively. 



and found a close-to-linear relation between a and C. 
Since the hydraulic resistance is -Rhyd = '^^hyd' 
elude that -R^yd depends linearly on CR^^^. Different 
classes of shape all display this linear relation, but the 
coefficients are non-universal. However, for each class 
only two points need to be calculated to fully specify 
the relation for the entire class. The difference is due 
to the smoothness of the boundaries. The elliptical and 
harmonic-perturbed classes have boundaries without any 
cusps whereas the rectangular and triangular classes have 
sharp corners. The over-all velocity profile tends to be 
convex and maximal near the center-of-mass of the chan- 
nel, see Fig. ^ If the boundary is smooth the veloc- 
ity in general goes to zero in a convex parabolic man- 
ner whereas a concave parabolic dependence is generally 
found if the boundary has a sharp corner (this can be 
proved explicitly for the equilateral triangle Q). Since 
the concave drop is associated with a region of low veloc- 
ity compared to the convex drop, geometries with sharp 
changes in the boundary tend to have a higher hydraulic 
resistance compared to smooth geometries with equiva- 
lent cross-sectional area. 

We believe that the explicit and simple link between 
i?hyd and C is an important observation since at the 
same time C is also central to the strength and effective- 
ness of various surface-related phenomena. We note that 
in micro-channels the flow properties and electrokinetic 
phenomena may be somewhat connected and substan- 
tial deviations from classical Poiseuille flow have been 
reported recently, see Ref. and references therein. Nev- 
ertheless, our observation is an important first step with 
relevance to the use of micro-fluidic channels in lab-on- 
a-chip applications. Furthermore, our results allow for 
an easy evaluation of the hydraulic resistance for ellip- 
tical, rectangular, and triangular cross-sections with the 
geometrical measure C being the only input parameter. 
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